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Abstract 

Applied to statistical physics models, the random cost algorithm enforces a Random 
Walk (RW) in energy (or possibly other thermodynamic quantities). The dynamics of this 
procedure is distinct from fixed weight updates. The probability for a configuration to be 
sampled depends on a number of unusual quantities, which are explained in this paper. This 
has been overlooked in recent literature, where the method is advertised for the calculation of 
canonical expectation values. We illustrate these points for the 2d Ising model. In addition, 
we prove a previously conjectured equation which relates microcanonical expectation values 
to the spectral density. 
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The performance of a numerical simulation depends on the chosen weight factors. Hence, 
some attempt should be made to optimize them for the problem at hand. The weight function 
of canonical MC simulations is exp{—f3E), where E is the energy of the configuration to be 
updated and f3 is the inverse temperature in natural units. The Metropolis algorithm and 
other methods generate canonical configurations {i.e. the Gibbs ensemble) through a Markov 
process. It has been expert wisdom [|l| for quite a while and became widely recognized in 
recent years that MC simulations with a-priori unknown weight factors are also feasible and 
deserve to be considered, for a concise recent review see 0]. For instance, weighting (in 
a certain energy range) with the inverse spectral density l/n{E) has turned out to be of 
practical importance. Examples are calculations of interfacial tensions for first order phase 
transitions, where improvements of many orders of magnitude were obtained. Instead of the 
energy other thermodynamic variables can be considered as well, e.g. 0]. To be definite, we 
focus on the energy. 

MC simulations with a-priori unknown weight factors require an additional step, not 
encountered in canonical simulations: A working estimate of the weight factors needs to be 
obtained first. Quite efficient recursive methods have been developed for this purpose 
Still, the question suggests itself whether one can possibly bypass the first step and develop 
methods which sample broad energy distributions right away. Indeed, it is possible to design 
updates such that a RW in some cost function is generated and the energy of a statistical 
physics model can be chosen as cost function. Unfortunately (as already noted in 0) the 
connection to the desired canonical expectation values is apparently lost. Nevertheless, it 
appears to be worthwhile to investigate properties of the thus generated configurations. In 
particular, some details are subtle and, besides the origin of the method, ignored in recent 
literature 0. 
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We consider generalized Ising models in d dimensions, described by the energy function 

E = - J2-JijSiSj (1) 

ij 

where the sum is over nearest neighbors and the exchange coupling constants Jij as well as 
the spins Sj, Sj take the values ±1. The Ising ferromagnet is obtained with Jij = 1. Other 
special cases are the Ising anti- ferromagnet, frustrated Ising models and spin glasses. We 
consider a configuration of N spins and choose periodic boundary conditions. Under flip of 
a single spin the energy can change by the following increments 

AEi = 4:i with i = -d, . . . , -1,0, 1, . . . ,d . (2) 

In the following we use i to label Flip Groups (FGs) of spins. We define now numbers 

d 

Ni with i = -d, -d + 1, ... d-1, d and ^ Ni = N (3) 

i=—d 

to partition the configuration of N spins with respect to the FGs. Namely, Ni denotes the 

number of spins which, when flipped, change the energy by AE = 4i. In the following A^, 
is referred to as Flip Group Magnitude (FGM). The random cost algorithm |^ achieves a 
RW in energy by flipping spins with suitable probabilities related to FGs. Let Pi be the 
probability for picking and flipping a spin in the FG labeled by A RW in E is obtained 
whenever the equation 

j2^{P^-P-^)=0 (4) 
i=l 

holds, because the expectation value of energy changes AE becomes then zero. In should be 
noted that Pq does not enter this equation and can be chosen at will. Besides, equation @ 
does not fix the other probabilities either, but allows considerable freedom concerning their 
further design. Before we come to this, let us note that has to be greater than zero for 
at least one i > 1 and one i < —1. Otherwise, a RW can no longer be achieved. This 
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latter difficulty happens in a local minimum or maximum of the system and, by whatever 
additional rule, one or more spins have to be flipped before the RW simulation can continue. 
In the following we assume that the noted A'^j > condition is fulfilled. 

Solutions to (I) are easily found, the following is given in p. We define 

A^+ = -^^iV,A^, and AE^" = -i- ^ iV„, A^_, (5) 
1=1 1=1 

where 

d d 

N+ = Y.N, and = ^ ^— (6) 

i=l i=l 

In the same way, we define 

d d 

P^ = Y.P^ and P- = Y,P^, . (7) 

i=l 1=1 

I.e. P"*" is the probability to pick any of the spins from the i > 1 FGs and P~ is the 
probability to pick any of the spins from the i < —1 FGs. Finally, assume that within those 
FGs the spins are picked with uniform probability, with = P^/N^ for i > 1 and with 
p~ = P~ /N" for i < —1. The RW equation (^) is then implied by the condition 

- P- AE- = P+ AE+ . (8) 

Choosing an arbitrary probability Pq, the probabilities P"*" and P~ follow immediately from 
this equation and the normalization condition Po + P^ + P~ = 1. Another way [^] to 
implement is to choose a spin at random and to reject the flip with the appropriate 
probability, then counting the configuration at hand again. Here we stay with (^. 

It follows from equations (^) and (§) that every such algorithm samples with weights 
which depend on the FGMs 

w = w {N_d, N_d+i, ... , Nd-i, Nd) . (9) 
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For configurations at a fixed energy E the implication of tliis equation is tliat the RW 
algorithm weights them differently depending on the FGM partition, whereas canonically all 
these configurations have the same weight. In the following we illustrate this point for the 
2d Ising ferro magnet. 

We have performed canonical as well as RW simulations for 2d Ising models on N = 
lattices with periodic boundary conditions. For the RW updating we used Pq = 0.2 and did a 
random flip, once the energy minimum was reached. At large energy we imposed a cut-off 
at E = 0, by replacing RW with random (canonical (3 = 1/(A;bT) = 0) updating for E > 0. 
To avoid getting lost in a flood of data, we focus on a single energy. After gaining some 
experience E/N = —1 with canonical simulations at /3 = 0.38 turned out to be a reasonable 
choice {E{P = 0.38) /N ^ —1). This value is in the disordered phase for /? below the Curie 
point at Pc = 0.51n(l + ^/2) = 0.4406 .... This correspond to a conflguration space region 
far away from the energy minimum E/N = —2 or the upper energy bound E/N = imposed 
on the RW simulation. The lattice sizes used are L = 4, 10, 20, 40 and 80. For each case we 
generated a statistics of 20 x 100, 000 sweeps through the lattice and calculated error bars 
with respect to twenty bins. For L = 4 we also obtained exact results by simply counting 
through all 2^^ conflgurations and convinced ourselves that the canonical simulation agrees 
(within very small statistical errors) with these exact results, whereas the RW simulation 
shows already considerable deviations. 

Let us focus on the microcanonical average values Ni/N. For the L = 80 lattice table 1 
collects results from the canonical as well as from the RW simulation. Although over-all 
small, in case of A^_2/A^ the discrepancy is about a factor of two and for all FGs the difference 
between the canonical and the RW values clearly exceed the error bars. In the average most 
spins, about 39%, are found in the FG with number i = 2. Hence, we have the best statistics 
for this FG and choose it to demonstrate a few more details. Although the discrepancy 
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between the N2/N values of the table appears quite small, there are considerable differences 
when we look at the distribution of N2. To correct for the expected (non-critical) finite size 
behavior, we define the quantity 

n = L-^ (N^ - N2) . (10) 
Figure 1 shows the canonical and the RW histograms h2{n), all normalized to 

L-'Y.h2{n) = l . 

n 

Error bars are negligible on the scale of this figure. The canonical histograms for the dif- 
ferent lattice sizes collapse nicely into one curve, whereas the finite size behavior of the RW 
histograms fails to reproduce this behavior. The RW histograms are far too broad and the 
peak height decreases with lattice size, implying that the discrepancy to canonical simulation 
increases with lattice size. The L = 4 data do not fit into the scale of this figure. As there 
are only four non-zero entries, we collect them in table 2. 

By taking averages, the dependence on the FGMs appears to be washed out. This is 
obvious for the average Ni/N values reported in table 1 and claimed to be true for 
thermodynamic quantities like the energy and the specific heat. The latter quantities were 
calculated from the spectral density g{E) which, in turn, was determined from the conjec- 
tured equations 

NiiE) g{E) = N^,{E + AE,) g{E + AE^) . (11) 

In appendix A we give a mathematical proof of these equations. From this it follows that 
the Ni{E) are microcanonical averages and not averages accumulated during the RW, as 
stated in 0. Only canonical, microcanonical or other simulations which give equal weights 
to distinct configurations at the same energy will converge towards the correct Ni{E) values 
for this equation. A rigorous calculation of the spectral density g{E) from RW data is not 
possible, because the RW weights depend on the FGMs (|^) 0. 
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For the L = 80 lattice table 1 lists RW and canonical expectation values for the Ni{E) 
at E/N = —1. As the differences are not too large, it is plausible that some RW and 
canonical results can be in qualitative agreement with one another. However, it is obvious 
that the dependence of the RW configurations weights (H) on the FGM partition {Ni} enter 
the Markov process. Even small deviations from the canonical weights may amplify, because 
they enter multiplicatively through each transition step. If the RW method is nevertheless 
used to estimate canonical expectation values, uncontrolled errors result with no guarantee 
that they will be negligible when it really matters (Murphy's law). 

A simulation is normally already subject to finite size and other difficult to control ap- 
proximations. Certainly, one would not like to build a large scale numerical investigation on 
a method which introduces an additional bias. The question arises, whether the weight de- 
pendence could eventually be controlled rigorously. Due to the large number of partitions 
of the total number of spins into FGMs (|^) the prospects for this do not look particularly 
good, but it may be worthwhile trying. Finally, we like to emphasize that the RW approach 
remains a competitive method for the purpose it was originally designated for, namely to 
find good energy minima for optimization problems and systems with conflicting constraints. 



Appendix A 

Here we prove equation (0) for the d- dimensional generalized Ising mode on a lattice with 
periodic boundary condition. Let K denote a spin configuration at energy E. By definition 

9[E) K 

holds. We introduce E' = E + AEi and label spin configurations at energy E' by K'. Then 
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and equation ([Tl|) becomes equivalent to 

EiV.W = E^-m- (12) 

K K' 

Equation (|12]) is shown as follows. We label spins in a fixed configuration by Sn{K) where 
n = 1, . . . ,N and is the (fixed) total number of spins. In this way each single spin is 
uniquely identified. The same is true for spins Sm{K') in the configurations at energy E' . 
Assume now, spin Sn^{Ki) is in flip group i and we flip it. It becomes a spin Smx{K'i) in the 
flip group —i at energy E' . No other spin, say s^jl-^a) with at least n2 ^ rii or Ki ^ K2, will 
be mapped on Smi{K[). The reason is: we can flip the spin Smi{K[) back and it will map 
precisely onto the orginal configuration and spin, i.e. become s„^(-ft'i). The same argument 
applies when we flip an arbitrary spin from flip group —i at energy E'. Together this proves: 
spins in the flip groups with magnitudes Ni{K) and N^i{K') are in on-to-one correspondence 
and, hence, equation ( |T^ ) is true. 
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Tables and Figure Captions 



i 


-2 


-1 





1 


2 


CS 


0.018853 (03) 


0.072752 (04) 


0.187630 (07) 


0.331070 (11) 


0.389694 (09) 


RW 


0.034282 (26) 


0.057936 (19) 


0.169412 (57) 


0.350240 (43) 


0.388130 (55) 



Table 1: Average Flip Group Magnitudes Ni/N, i = —2, ... , 2 as obtained for the Canonical 
Simulation (CS) versus the Random Walk (RW) simulation on a 80 x 80 lattice. Error bars 
are given in the parenthesis and apply to the last two digits. 



Figure 1: Canonical and RW histograms h2{n) at E/N = —1 with n given by (plQl). All 
the canonical histograms collapse onto the highest curve. The RW histograms follow 
in the order L = 10, 20, 40 and 80 from up to down. 
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n 


-1.5276 


-0.5276 


-0.0276 


0.4724 


N2 





4 


6 


8 


Exact 


32/424 = 


192/424 = 


1088/424 = 


384/424 = 


Exact 


0.07547... 


0.45283... 


2.56603... 


0.90566... 


CS 


0.0759(11) 


0.4543(35) 


2.5672(43) 


0.9026(31) 




0.221(3(08) 


0.5171(07) 


2.3080(11) 


0.9233(10) 



Table 2: The h2{n) histogram results for the L = 4 lattice (on this lattice there are 424 
configurations at E/N — —1). 
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